arXiv: 1507.04957v2 [cond-mat.str-el] 30 Sep 2015 


Hard-core boson approach to the spin-1/2 triangular-lattice antiferromagnet Cs 2 CuCl 4 
at finite temperatures in magnetic fields higher than the saturation field 
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We study the high magnetic field regime of the antiferromagnetic insulator CS2CUCI4 by ex¬ 
pressing the spin-1/2 operators in the relevant Heisenberg model in terms of hard-core bosons and 
implementing the hard-core constraint via an infinite on-site interaction. We focus on the case 
where the external magnetic field exceeds the saturation field B c « 8.5 T and is oriented along the 
crystallographic a axis perpendicular to the lattice plane. Because in this case the excited states are 
separated by an energy gap from the ground state, we may use the self-consistent ladder approxi¬ 
mation to take the strong correlations due to the hard-core constraint into account. In CS2CUCI4 
there are additional interactions besides the hard-core interaction which we treat in self-consistent 
Hartree-Fock approximation. We calculate the spectral function of the hard-core bosons from which 
we obtain the in-plane components of the dynamic structure factor, the magnetic susceptibility, and 
the specific heat. Our results for the specific heat are in good agreement with the available ex¬ 
perimental data. We conclude that the self-consistent ladder approximation in combination with a 
self-consistent Hartree-Fock decoupling of the non-hard-core interactions gives an accurate descrip¬ 
tion of the physical properties of gapped hard-core bosons in two dimensions at finite temperatures. 

PACS numbers: 75.10.Jm, 05.30.Jp, 75.40.Cx, 75.40.Gb 


I. INTRODUCTION 

The magnetic behavior of the antiferromagnetic insu¬ 
lator CS2CUCI4 can be described by a spin-1/2 Heisen¬ 
berg model on a triangular lattice where the interlayer 
coupling is much weaker than the intralayer couplings. 
Due to the relatively weak exchange couplings, a field 
induced ferromagnetic ground-state can be reached at 
fields larger than the saturation field B c ss 8.5 T 
where the magnetic field is along the crystallographic a 
axis perpendicular to the lattice plane. This allows a 
precise measurement of the exchange couplings via in¬ 
elastic neutron scattering experiments where the single 
magnon dispersion gives direct access to the exchange 
couplings^ CS2CUCI4 has been intensively studied due 
to its interesting properties, e.g., spin-liquid behavior 
with spinon excitations }- 2 ® ^ Bose-Einstein condensation 
of magnons at the quantum critical point EEMUl an d a 
rich phase diagram for in-plane magnetic fieldsPESlIMlA. 
diverse range of observables have been experimentally 
investigated: dynamic structure factor,®^ electron spin 
resonance spectra}^! magnetic susceptibility}^ magneto¬ 
caloric effect}^ 6 nucl ear magnetic resonance relaxation 
ratepl specific heatp^ and ultrasound velocity and 
attenuation. 27 ^^ 

In this work, we consider the case of a large magnetic 
field B > B c along the a axis, where the magnon excita¬ 
tions are gapped and the ground state is the fully mag¬ 
netized ferromagnet. Our goal is to describe the thermal 
excitations above the ground state and to compare with 
experimental results for the specific heatP^GI] We base 
our theoretical approach on a mapping of the spin-1/2 
operators to hard-core bosons EMU p or magnetic fields 
B > B c and low temperatures, we then have a dilute gas 


of gapped hard-core bosons where the ladder approxi¬ 
mation captures the leading order low-temperature con¬ 
tributions to the self-energy. Although the ladder ap¬ 
proximation has been extensively applied to the Bose- 
condensed phase of dilute gases of hard-core bosons (see 
Ref. [T5] and references therein), some subtleties related 
to the hard-core limit in the gapped phase have only re¬ 
cently been discussed by Fauseweh, Stolze, and Uhrig 
(FSU)P 2 33 Benchmarking the ladder approximation for 
an exactly solvable one-dimensional model of hard-core 
bosons, FSU found that the ladder approximation in¬ 
deed reproduces the correct low-temperature behavior 
and that a self-consistent ladder approximation even ex¬ 
tends the applicability to arbitrarily high temperatures.^ 
In this work, we apply the self-consistent ladder ap¬ 
proximation to the relevant two-dimensional model for 
CS2C11CI4. For a realistic description of this material we 
have to include additional interactions apart from the 
infinite on-site interaction describing the hard-core con¬ 
straint. To further explore the range of validity of the 
self-consistent ladder approximation, we have also ap¬ 
plied this method to the exactly solvable one-dimensional 
XY model; extending the analysis of FSUj22 we have ex¬ 
amined the breakdown of the self-consistent ladder ap¬ 
proximation in the vicinity of the quantum critical point 
of this model. 

The rest of this work is organized as follows. In the 
next section, we introduce the relevant spin model for 
Cs 2 CuCl4 and describe the mapping of this model to 
an effective hard-core boson model. Then, in Sec. |III[ 
we describe our theoretical approach based on the self- 
consistent ladder approximation for the hard-core inter¬ 
action and a self-consistent Hartree-Fock decoupling for 
the remaining non-hard-core interactions. In Sec. m we 
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FIG. 1. (Color online) Part of the anisotropic triangular lat¬ 
tice formed by the spins of CS 2 CUCI 4 . The stronger exchange 
coupling J connects nearest-neighbor spins along the crys¬ 
tallographic b axis, while the weaker exchange coupling J' 
connects nearest-neighbor spins along the diagonals. There 
are also weak Dzyaloshinskii-Moriya interactions D = dcDz 
connecting neighboring spins along the diagonals where the 
direction of D is indicated by © for +z and © for —z. We 
consider only the case where the magnetic held B = Bz is 
along the a axis perpendicular to the plane of the lattice. 


investigate the breakdown of the ladder approximation 
near the quantum critical point for the exactly solvable 
one-dimensional XY model and in Sec. [V] we present 
our numerical results for CS 2 CUCI 4 , which we compare 
with experimental data for the specific heat. Finally, in 
Sec. |VI| we summarize our main results. In three appen¬ 
dices we give additional technical details of our calcula¬ 
tions. 


II. HARD-CORE BOSON MODEL FOR 

Cs 2 CuCl 4 


Here, x , y, and z are the unit vectors of our Carte¬ 
sian coordinate system. Due to the fact that inver¬ 
sion symmetry is broken for CS 2 C 11 CI 4 , there are also 
Dzyaloshinskii-Moriya (DM) interactions D l:) = D(Ri — 
Rj)z connecting neighboring spins along the diagonals 
with Z?(±£ 2 ) = £>(±<$ 3 ) = TD. The precise form of 
the Hamiltonian (2.1 1 and the values of the interaction 


constants have been measured by inelastic neutron scat¬ 
tering experiments in magnetic fields higher than the sat¬ 
uration field B c = 8.44(1) T. The accepted values areP 
J = 0.374(5) meV = 4.34(6) K, J'/J = 0.34(3), and 
D/J = 0.053(5). There is also a weak interlayer cou¬ 
pling J"/J = 0.045(5), which we neglect because it is 
only important at very low temperatures T < 0.1 K and 
in the antiferromagnetically ordered phase in magnetic 
fields B < B c which we do not consider here. Recently, 
additional DM interactions, including in-plane compo¬ 
nents, have been measured via electron spin resonance 
experiments.^ We neglect these additional DM inter¬ 
actions because they are mainly important for in-plane 
magnetic fields.® Furthermore, our theoretical approach 
relies on the U( 1) symmetry due to the spin-rotational 
invariance with respect to the z axis, which would be 
broken by in-plane DM interactions. 

In this work, we will use the har d-co re boson repre¬ 
sentation of the spin-1/2 operators!^*® Recall that the 
spin-1/2 operators fulfill the commutation relations 

[St,Sj] = 2 SijSt, [Sf , S*] = zpSijSf , (2.3) 

where Sf = Sf ± iSf and Sf = 3/4. Additionally, the 
spin-1/2 operators obey an on-site exclusion principle}® 


S+S- + S~S+ = 1 , (S+y = ( S-y = 0 . (2.4) 


It has been established that the magnetic behav¬ 
ior of CS 2 CUCI 4 can be described by the following 
two-dimensional antiferromagnetic spin-1/2 Heisenberg 
model in an external magnetic field along the crystallo¬ 
graphic a axis,® 


To realize these relations, we can express the spin opera¬ 
tors in terms of hard-core boson creation and annihilation 
operators, 

S+ = bi, S~=bl Sr = 1/2-b\bi, (2.5) 


n - \ Y. J " s ‘' + D ’ ■ uV ' x si)l - 

ii i 

( 2 . 1 ) 

where the summations run over all N lattice sites, h = 
g^BB is the Zeeman energy associated with an external 
magnetic field B = Bz, and g = 2.19(1) is the effective 
g- factor.3The spin-1/2 operators Si = S(Ri) are located 
at the lattice sites Ri of an anisotropic triangular lattice 
with lattice constants b and c, as shown in Fig. [T] The 
exchange couplings Jij = J(Ri — Rj) connect nearest 
neighbors along the crystallographic b axis and along the 
diagonals with J(±<$i) = J and J(±<$ 2 ) = J(±S 3 ) = J', 
where the three elementary lattice vectors are 


where the hard-core boson operators satisfy the commu¬ 
tation relation 


h. h 1 
o u o 3 


= Sij (l - 2 b\bi) , 


( 2 . 6 ) 


and the occupation number per site is restricted to hi = 0 
or 1, with hi = b\bi. The hard-core boson constraint 
and the commutation relation ( 2 . 6 ) can be realized by 


treating the hard-core bosons as canonical bosons with 
an infinite on-site repulsion, 


U 


Ru = 7 T b\b\bibi, with U -A 00 . (2.7) 


<5i = bx, S 2 


b „ c ^ 
-2 X+ 2 V ’ S ‘ 



(2.2) 


Note that the magnon excitations of the underlying spin 
system correspond to hard-core boson excitations. 
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Using Eq. (2.5) to express the spin operators in our 
Hamiltonian (2.1) in terms of hard-core bosons, we obtain 
the following hard-core boson Hamiltonian, 


6 — + 2jy (Jq + U)b\ l+q b\,,_ q bk , bk +Eq, 

k k,k' ,q 

(2.8) 

where we have Fourier transformed the hard-core boson 
creation and annihilation operators, 


bk 


1 

Vn 


J2 bie~ ik R ’ 



1 

Vn 


J2 b 


t ik-Ri 


i 


(2.9) 



In the following, we will neglect the unimportant con¬ 
stant energy term 


= N (f - i) ■ < 2 - 10 > 

The excitation energy in the quadratic part of the 
Hamiltonian can be written as 

6c = £fc - /b (2-11) 

where we have introduced the chemical potential 

[i = h c — h, (2.12) 


FIG. 2. (Color online) Contour plot of the energy dispersion 
£fc defined in Eq. (2.131. The white cross marks the minimum 
of £ k at Q « (3.4747M). 


the experimental value is more accurate than a calcula¬ 
tion via the Hamiltonian (2.8 1 . The reason is that the 
interaction constants in Eq. (2.8 1 have some experimen¬ 
tal uncertainty and we have also neglected the interlayer 
coupling J"\ including J" in the calculation would re¬ 
sult in B c ss 8.5 Tj 15 * 1 '* The value of the saturation field 
B c is important because, for a given magnetic field, it 
determines the energy gap 


and the energy dispersion 


D\ 


£ k ~ 2 (Jk ~ Jq 


Here 


A = —jj, = h — h c . (2.18) 

(2.13) We note that a small change of B c by 0.04 T changes the 
gap by about 0.014 J, which is only significant close to 
the quantum critical point. 


Jk — 6c iDki 


(2.14) 


where the Fourier transforms of the exchange and III. IMPLEMENTING THE SELF-CONSISTENT 
Dzyaloshinskii-Moriya interactions are LADDER APPROXIMATION 


■4 


Dk 


Y^J{R)e~ ikR 

R 

2 J cos (k x b) + 4 J’ cos 




J2D(R)e~ ikR 

R 


—AiD sin 



cos 



(2.15) 


(2.16) 


In Eq. (|2T3l, J R 


—2.325 J is the absolute minimum 


of Jk at Q S3 (3.474/6,0). Finally, the saturation field is 
given by 


Be = — = ( J 0 D - J%) « 8.4 T. (2.17) 

g^B Zg^B 

A contour plot of £k is shown in Fig. [2j In the following, 
we will use the direct experimental value of the satura¬ 
tion field B c = 8.44(1) T instead of B c R3 8.4 T because 


In this section, we explain our theoretical approach 
to the hard-core boson Hamiltonian (2.8 1 . The central 
problem is how to deal with the interaction part of the 
Hamiltonian, 


6 int = ^E( J ? + U)bl +q bl_ q bvb k , ( 3 . 1 ) 

k,k' ,q 

containing the exchange interaction J q and the infinite 
hard-core interaction U —> oo. We will deal with both 
interactions using different methods: for the J q part we 
use a self-consistent Hartree-Fock decoupling, while for 
the hard-core interaction U we use the self-consistent lad¬ 
der approximation P3 This is necessary because the J q in¬ 
teraction cannot be easily included in the self-consistent 
ladder approximation, as this would not allow a direct 
solution for the effective interaction F from the Bethe- 
Salpeter equation which would significantly complicate 
matters, especially regarding the limit U —> oo. 




























4 


A. Hartree-Fock decoupling 


We approximate the J q interaction term in Eq. (3.1) 
using a self-consistent Hartree-Fock decoupling. There¬ 
fore, we write this term in real space, 


^ E JA + A'- q b ^ = \T, J ^ b Mb j , (3.2) 

k.k' ,q i,j 

and then we apply the usual Hartree-Fock decoupling, 
blbibjbj ss riib^bj + rijb\bi — riiUj 

+ rubify + Tijbjbi - TijTji, (3.3) 

giving 

2 E Ji j h \ b i h ) h 3 ~ E (".M 1 ’' + T H b \ b i) + £ MF, 

ij i,j 

(3.4) 

where the Hartree-Fock parameters are given by 

rii = (b\b^ , Tij = (b\bj^j , (3.5) 

and the constant energy term is 

Emf = - 9 E Ji J ( n * n f + T H T ji) ■ (3-6) 

Due to translational invariance, we have 

Tli 1 7., Tij T(R t Rj). (3.7) 

Because there is no inversion symmetry, t(R) is a com¬ 
plex number satisfying t*(R ) = t(—R). Since the ex¬ 
change coupling J(R) is only non-zero for R = ±<5,, 
there are three complex Hartree-Fock parameters related 
to t(R), 


ti=t(S i), r 2 = r(5 2 )j t 3 = t(5 3 ). (3.8) 


However, the Hamiltonian (2.8 1 is invariant under the 
transformation k y —> — k y and therefore r 2 = t 3 . Trans¬ 
forming Eq. (3.4) back to momentum space, we get 


9jy E Jqhk+q^k'-q b k' b k ~ E] + n ^o) ^fc^fc+^MF, 

k,k' ,q k 

(3.9) 

where 


temperature; moreover, the Hartree-Fock approximation 
leads to a renormalization of single-particle excitation en¬ 
ergies £& —> where the renormalized excitation ener¬ 

gies are 

ik — £k — A 4 + + nJo- (3.12) 

The self-consistency equations for the Hartree-Fock pa¬ 
rameters are given by 


n = r(Si) = ^^n k e tk ' Si , (3.13a) 

k 

n=-^^n fc , (3.13b) 

k 

where the occupation number of a state with momentum 
k is given by 

n k = (b{b k ^ . (3.14) 


If we neglect the hard-core interaction, we simply obtain 
the Bose-Einstein distribution, 


nT = 


1 


e /3ff fc _ 1 ’ 


(3.15) 


where /3 = 1/T is the inverse temperature and the 
renormalized excitation energy £*. can be obtained in a 
straightforward way by solving the self-consistency equa¬ 
tions for the Hartree-Fock parameters for U = 0. Ne¬ 
glecting the hard-core interaction is possible only for 
small temperatures T <C J when the bosons are so di¬ 
lute that the hard-core interaction does not contribute 
significantly. 


B. Self-consistent ladder approximation 

After the Hartree-Fock decoupling of the J q interac¬ 
tion, we obtain a Hamiltonian where the only remaining 
interaction is the infinite on-site repulsion, 

n = E Zk b l b k + E b i+q b i' -qbk'bk + E MF . 

k k,k' .q 

(3.16) 

We will deal with this hard-core interaction using the self- 
consistent ladder approximation developed in Ref. [52] . 


1. Imaginary time path integral formalism 


Jl = 2JRe (r ie ifc ' 51 ) + 2J'Re (r 2 e ,;fc <S2 + r 3 e ikS3 ) 


Emf = —N 


■y n 2 + ^l r il“ + J’ (l r 2|“ + |d| 2 ) 


(3.10) 


. (3.11) 


To derive the self-consistent ladder approximation, it 
is convenient to formulate the problem in terms of an 
imaginary time path integral.^ The E uclidean action as¬ 
sociated with the Hamiltonian (3.161 is 


S[b,b] = - Gq 1 (K)b K b K 


ik 


U 


bK+Q b K'-QbK' b K■ (3-17) 


The Hartree-Fock approximation gives a constant energy 
shift F/mF) which depends on the magnetic field and the 


K,K',Q 
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FIG. 3. Diagrammatic representation of the se lf-ene rgy in 
terms of the effective interaction as given in Eq. (3.271. 


FIG. 4. Ladder approximation for the effective interaction T 
including all particle-particle ladder diagrams. 


Here, we have introduced the composite index I\ = 
(k. iuJk) with the corresponding sum 





(3.18) 


where uj k are bosonic Matsubara frequencies. The com¬ 
plex boson fields in imaginary time have been Fourier 
transformed to frequency space as 

W = -^5>- iWfeT fc*, (3.19a) 

- bk{T ) = -±=Y,e iUkT ~ b K- ( 3 - 19b ) 

ivfc 

The Green function G(K) and the corresponding self¬ 
energy E(/\) are defined via the functional average 


(b K b K ) = -pNG(K) 


■/3N -,-, 

Gq 1 (K) — E(A') 


(3.20) 


where the bare Green function Gq(K) is given by 


G 0 (K) = , 1 . (3.21) 

— ffc 

From this path integral formalism a perturbative dia¬ 
grammatic expansion of the one-particle irreducible self¬ 
energy E(/\) can be obtained in terms of the bare Green 
function Gq(K) and the interaction U. 


2. Self-consistent ladder approximation 



FIG. 5. Diagra mmatic representation of the Bethe-Salpeter 
equation (3.241 for the effective interaction T. 


ladder diagrams only include a single line of this type; 
each such line gives a suppression of exp(—/3A). This 
can be seen by considering the bare Green function in 
imaginary time, 


G 0 (fc,r) 



(3.22) 


where ns{x) denotes the Bose function, 

= ( 3 - 23 ) 


We see that ns(^fc) oc exp(—/3A) for /3A 1 and there¬ 

fore Go(k, t) a exp(—/3A) for r < 0. 

We can go beyond the ladder approximation and in¬ 
clude higher order terms by using the full Green function 
G(K) in the diagrammatic expansion instead of the bare 
Green function Gq(K) and then finding a self-consistent 
solution. In this self-consistent ladder approximation, 
the effective interaction fulfills the Bethe-Salpeter equa¬ 
tion shown in Fig. [5] 


Since we are dealing with a strictly non-perturbative 
problem (U —i oo), it is necessary to sum over a suitable 
infinite set of diagrams containing infinite powers of U. 
Here we approximate the self-energy by summing over 
all particle-particle ladder diagrams, where we express 
the self-energy in terms of the effective interaction T, as 
shown in Fig. [3] The effective interaction then includes 
the infinite series of particle-particle ladder diagrams in¬ 
dicated in Fig. [4] Formally, this approximation is jus¬ 
tified for /3A 3> 1 because the neglected diagrams are of 
order exp(—/3A) smaller than the ladder diagrams. The 
neglected diagrams for the self-energy include at least 
two lines going backwards in imaginary time while the 


r(A", K-P) =U-U F(A'", K\ P)G(K")G(P — K"). 
J K" 

(3.24) 

Because the hard-core interaction U is a constant inde¬ 
pendent of the momentum transfer, the Bethe-Salpeter 
equation has the simple solution 

T(K',K; P) = r(p) = Y+JJn(pj’ <^> 

where we have defined the particle-particle bubble 

n(P) = [ G(Q)G(P — Q). (3.26) 

Jq 
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The self-energy in the self-consistent ladder approxima¬ 
tion is given by 


£(if) = -2 / G(Q)T(Q + K)e l “ q0+ , (3.27) 

JQ 

which is shown diagrammatically in Fig. [3] The conver¬ 
gence factor e lu>q0+ implements the correct time ordering 
at the interaction vertex in the first order term in U where 
a propagator line starts and ends at the same vertex!-® 
The Green function has the spectral representation 

G(K) = r dx A ^ k ’ x \ ( 3 . 28 ) 

J- oo - x 


where the spectral function is given by 


A(k, oj) = -ImG(fc, oj + i 0 + ) 


7T 

1 


Im£' R (fc, oj) 


oj — — ReY, R (k,oj) + [Im£ fl (fc, w)]“ 

(3.29) 

and the retarded self-energy is obtained by analytic con¬ 
tinuation to real frequencies, 


£ fl (fc,w) = £(fc,w-H0+). 


(3.30) 


We note that the spectral function of hard-core bosons 
fulfills the following sum rulef® 


dojA(k,oj) = 


b k ,b{ 


= 1-2 n. 


(3.31) 


with 

1 

N 

Our goal is to calculate the self-consistent solution for 
the spectral function A(k, oj). But before we can do this, 
we have to take the limit U —> oo analytically. 



/ OO 

dx A{k, x)n B {x). (3.32) 

-oo 


3. Taking the limit U —> oo 


We can write the effective interaction T(P) as 


T(P) = 

where the second term 
ST(P) = 


n(P) 


WT(P), 


n(P) l + pn(P) 


(3.33) 

(3.34) 


does also contribute to the U —i oo limit because the 
denominator of ST(P) can vanish at high frequencies 
oj p ~ 0(U) leading to an additional delta function con¬ 
tribution which has to be taken into account. This sub¬ 
tlety of the U —> oo limit has been noticed only quite 


recently by FSU.^ We now follow FSU to derive the cor¬ 
rect hard-core limit for our model. First of all, we note 
that II(p, oj) oc 1/oj for oj —> oo. This allows us to intro¬ 
duce the spectral representation 

n (p,oj) = r dx (3.35) 

J- oo W - X 


where 


p(p, oj) =-ImII(p, oj + z0 + ) 

7T 

1 f°° 

= 5Z / dxA(q,x)A(p-q,oj-x) 

iV q d-oo 

(3.36) 


x [n B {x) - n B (-x)\. 

Now, we use the fact that for oj —> oo 

rl ^ tt ~ u2u (P^) 1 

T p,w) - U = -r oc -, 

l + UU(p,oj) oj 


(3.37) 


because II(p, oj) oc 1/w for u -> oo. This implies that 
r(p, oj) — U has the spectral representation 


r( P ,u)-u = 


oj - x 


(3.38) 


where 

P(P,u) = ^“ Im [ r ( 2 h w + * 0+ )] 

= f(P , w ) - hm [< 5 r ( p , oj + z 0 + )] , ( 3 . 39 ) 

with 


f(p,u) = - 


~p(p,oj) 


V f°° dx 

J — OO 


p(p,x) 


[ttp(p,oj)Y 


(3.40) 


Here, V denotes the Cauchy principal value which arises 
from the identity l/(w-|-z0 + ) = V(l/oj) — iirS(oj). For the 
contribution of <5T to p(p, oj) we recall that the denomi¬ 
nator of <5r(p, oj) can vanish when oj ~ 0(U) and only in 
that case there can be a contribution from 5T. Therefore, 
we expand n(p, w) for large frequencies oj ~ 0{U) (we 
take U to be very large but finite), 


H {p,oj) 


Poip) + Pi(p) { c 

OJ OJ 2 



(3.41) 


where 


/ OO 

dxp{p,x), (3.42a) 

-OO 

/ OO 

dxxp(p,x). (3.42b) 

-OO 


We find 


6T(P) « - 


Mp) + PiiPl +0 /i\ 

UJ 1 UJ 2 XUJ' 5 / 

1 


1 + jjPgM + [//tUPJ + Of 4) 

UJ UJ 2 V UJ^ / 


pi(p) 


(3.43) 
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where the terms Oil/oo 3 ) and 0(U/co 3 ) vanish at the pole 
co ~ 0(11) in the limit U —> oo, justifying the expansion 
to order 1/co 2 . Therefore, we have 


ST(P) « - 


u> 


Mp) + pi(p) ( w - W i(p)) (uj - uj 2 {p )) ’ 

UJ UJ Z 

(3.44) 

where the poles are given by 


U 1 ( P ) = ~ UP0(P) 


2 V 4 
-Up 0 (p), U ^ oo 


- ^Pi(p) 


(3.45a) 


w 2 (p) = - 


^Po(p) | JU 2 Pq(p ) 

Pi(p) 


- ^Pi(p) 


Po(p)’ 


U —► oo. 


(3.45b) 


Only the pole at uq ~ 0([/) is relevant for the ana¬ 
lytic continuation in Eq. (3.391 because the other pole 
at u )2 ~ O(U 0 ) is spurious, since we have expanded for 
large frequencies co ~ 0(U). In total, we get 


p(P,v) = f{P,v)~ 


1 


p°(p) + Mp) u - W2 (p) 


<5(w-wi(p)). 
(3.46) 


We can now use the spectral representation (3.381 in 


Eq. (3.271 and take the limit U —> oo to get the following 
expression for the self-energy, 


r\ /»00 f OO 

E(A') = -—dx dx'A(q, x')f(q + k, x) 

„ J— oo J — OO 


n B [x) - n B (x') 
ico k + x' — x 


4e 


dx A{q 1 x)n B (x) 


x + iujk pi(q + k) 


p Q {q + k) pl(q + k) 


(3.47) 


By analytic continuation to real frequencies, we obtain 
the real and imaginary part of the retarded self-energy, 

2 _ r°° 

Re£ fl (fc,w) = ^ ^ / dx A(q,x)n B (x) 

^ J ~ OO 


i + w pi(q + fe) 


_p 0 (<? + fc) Po(9 + fe). 

PT.{k,x) 


POO 

+V / dx 


CO — x 


ImE fl (fe,w) = —irpz{k,co), 


(3.48a) 

(3.48b) 


where 


2 _ r°° 

P-z(k,u) =—^2 dx A(q,x)f(q + k,x + uj) 

^ „ J — OO 


x [n B (x) -n B (x + w)]. 


(3.48c) 


To summarize, we have obtained the self-energy in the 
limit U —> oo which we can calculate starting from an 
initial spectral function. Via Eq. (3.291 we can then cal¬ 


culate the next iteration of the spectral function allowing 
us to find a self-consistent solution for the spectral func¬ 
tion. After each iteration the Hartree-Fock parameters n, 
Ti, t 2 , and r 3 have to be updated via the self-consistency 
equations (3.13a I and (3.13b) using 


n k = ( b\b k ) = 


dx A(k , x)n B (x). 


(3.49) 


IV. NUMERICAL RESULTS FOR THE 
ONE-DIMENSIONAL XY MODEL 

Before applying the above approach to the hard-core 
boson model for CS 2 C 11 CI 4 , it is instructive to test its va¬ 
lidity for the exactly solvable one-dimensional spin- 1/2 
XY model in a magnetic field. Although a similar model 
has been already studied in detail in Ref. [32J, the break¬ 
down of the self-consistent ladder approximation in the 
vicinity of the quantum critical point has not been inves¬ 
tigated. 

The XY model in one dimension is given by 

Km = JJ2 (W+i + S i S i+i) - h E » (4- 1 ) 

i i 

which can again be mapped to hard-core bosons (neglect¬ 
ing constant terms), 


= J2 ^b\b k + E Y, b l + A'- q h' b k, (4.2) 

k k,k' ,q 

with excitation energy 

£fc = J [cos (k x b) + 1 ] - p, (4.3) 

where p = h c — h = — A and h c = J. An exact solu¬ 
tion can be found by mapping the hard-core bosons to 
fermions via the Jordan-Wigner transformation, 

bj = c ‘ Cl Cj, b] = C ‘ CI , (4.4) 

resulting in the quadratic Hamiltonian 

%id = E&4*. (4-5) 

k 

where the operators cj) and c k are fermionic creation and 
annihilation operators. The hard-core boson density is 
therefore exactly given by 

"“vEis/yc < 4 - 6 > 

k 

which can be compared to the approximate solution from 
the self-consistent ladder approximation for the hard-core 



































FIG. 6 . (Color online) Comparison of the results for the boson 
density n obtained from the self-consistent ladder approxima¬ 
tion (symbols) with the exact result (solid lines) at different 
energy gaps A for the one-dimensional XY model. 
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FIG. 7. (Color online) Relative error Sn/n of the boson den¬ 
sity (see Eq. (4.8 1 ) in the self-consistent ladder approximation 
for the one-dimensional XY model as a function of tempera¬ 
ture at different energy gaps A. 


boson model. We note that for low temperatures at /r = 
0, the exact density (4.6 1 has the following asymptotic 
behavior, 


dx^J— « 0.241 y/Tfj, (4.7) 

7T Jo e + 1 

in agreement with the expected behavior of one¬ 
dimensional bosons at the quantum critical point P® 

In Fig. [ 6 ] we compare the approximate result for the 
boson density n obtained within the self-consistent ladder 
approximation with the exact solution. The relative error 
Sn/n of the approximate result is shown in Fig. [7] where 
we define 


Sn/n 


^ladder ^exact 


^-exact 


(4.8) 


Here, redder is the result from the self-consistent ladder 
approximation and n ex act is the exact result. For a finite 
gap, we see that the error vanishes for low and high tem¬ 
peratures with a maximum error at T ss A, while the 


error becomes smaller for larger energy gaps. For A = 0, 
the error keeps increasing for smaller temperatures, get¬ 
ting closer to the quantum critical point at T = 0, but 
decreases for higher temperatures. We conclude that the 
self-consistent ladder approximation gives good results 
over the whole temperature range for finite energy gaps 
A >0.1 J. 

V. NUMERICAL RESULTS FOR Cs 2 CuCl 4 

In this section, we apply our hard-core boson approach 
described in Sec. |III| to the relevant model for CS 2 CUCI 4 
given in Sec. [TTJ From the numerical solution of the self- 
consistent ladder approximation we calculate the spec¬ 
tral function of the hard-core bosons at finite tempera¬ 
tures for different magnetic fields in the regime B > B c 
where the energy gap A > 0 is finite. Given the spec¬ 
tral function, we can calculate the magnetization, the 
internal energy, and the transverse part of the spin dy¬ 
namic structure factor. From the magnetization and in¬ 
ternal energy we obtain the magnetic susceptibility and 
the specific heat by numerical differentiation. Finally, we 
compare our results with experimental data for the spe¬ 
cific heat ! 1,r>|16 1 Technical details of the numerical solution 
of the self-consistent ladder approximation can be found 
in Appendix A and Appendix B. 


1. Spectral function 

Due to the finite energy gap A > 0, at zero temper¬ 
ature the spectral function is exactly given by the non¬ 
interacting spectral function, 

A{k,u) = A 0 (k,u) = 5(u) - £ fc ). (5.1) 

At finite temperatures, interactions will lead to a renor¬ 
malization of the excitation energy and a broadening 
of the delta peaks. Since we are treating the J q interac¬ 
tion term on a Hartree-Fock level, this alone would only 
renormalize the excitation energy by £& —> £& resulting 
in 


A(k, u) = S(lo — £&). (5.2) 

Taking in addition the hard-core interaction via the self- 
consistent ladder approximation into account will lead to 
a broadening of the spectral function with rising tem¬ 
perature, as shown in Fig. [ 8 ] Besides the broadening, 
we notice that the bandwidth shrinks with rising tem¬ 
perature and the minimum of the spectral function gets 
shifted to higher energies increasing the effective energy 
gap from its bare value A at T = 0. In Fig. [9] we con¬ 
trast the behavior of the spectral function at k = 0 and 
at the minimum of the dispersion k = Q. While at k = 0 
the position of the peak only moves to slightly higher en¬ 
ergies, at the minimum of the dispersion the peak gets 
considerably shifted to higher energies. Due to the finite 
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frequency resolution in our numerical calculation (see 
Appendix A), we cannot reach arbitrarily low tempera¬ 
tures and are restricted to temperatures T > 0.2A where 
the spectral function is not too narrow to be resolved. 
However, in the temperature range T < 0.2A, the hard¬ 
core interaction can be neglected and we can then just 
use the self-consistent Hartree-Fock decoupling without 
hard-core interaction, as we will show further below. 

The spectral function can be related to the in-plane 
components of the spin dynamic structure factor. The 
spin dynamic structure factor is defined by 


where n is the boson density per site which can be ex¬ 
pressed in terms of the spectral function, 

n = jjJ2 nk ’ ( 5 - 9 ) 

k 

with 

n k = (b\b k \= f dx A(k, x)n B {x). (5.10) 

' J —OO 

We define the magnetic susceptibility \ via 



(si k {t)Si{ 0 )), 


(5.3) 


where a, /3 = x,y,z and the Fourier transforms of the 
spin operators are defined via 


qa _ \ ' -ik-Ri qo 

_ y/N . 6 


(5.4) 


The in-plane components of the spin dynamic structure 
factor are given by 


S xx {k,oj) = 


S xy (k, oj) = 


(5.5b) 

where S xy (k,oj) and S yx (k,uj) do not vanish due to the 
broken inversion symmetry. The U( 1) symmetry due to 
the spin-rotational invariance with respect to the z axis 
requires that 

S xz (k,tu) = S zx (k,Lo) = S yz (k,co) = S zy {k,cu) = 0. 

(5- 6) 

The S zz component of the spin dynamic structure factor 
cannot be simply expressed in terms of the spectral func¬ 
tion because it is a two-particle Green function in terms 
of the hard-core boson operators, 


S vv (k, u) 

4 l _ e -pu W ) + A (~ k > W )) : 


(5.5a) 


- S yx (k,uj ) 
1 1 
41 - e~P u 


(A(k,u) - A(-k,w )), 


dm 

dB' 


(5.11) 


In the limit T —> 0 the boson density vanishes and the 
asymptotic low-temperature behavior of the susceptibil¬ 
ity is therefore the one of free bosons because all interac¬ 
tion are frozen out, 


X 


oc 


T ^ e -A/T 


(5.12) 


where d is the dimensionality {d = 2 in our case) and A 
the energy gap. 

The numerical results for magnetic moment and mag¬ 
netic susceptibility for different magnetic fields above the 
saturation field are shown in Fig. [TO] and Fig. m respec¬ 
tively. In Fig. |12[ we compare our numerical results from 
the self-consistent ladder approximation with the low- 
temperature Hartree-Fock approximation without hard¬ 
core interaction and with a simple spin mean-field theory, 
which we describe in Appendix C. We see that for low 
temperatures T <C J the Hartree-Fock and the ladder 
approximation give essentially the same results. This al¬ 
lows us to use the Hartree-Fock approximation in the 
low-temperature regime where the self-consistent ladder 
approximation is difficult to implement due to the lim¬ 
ited frequency resolution. At higher temperatures the 
hard-core interaction becomes important and the high- 
temperature behavior is approximately captured by the 
spin mean-field theory. 


3. Internal energy and specific heat 
The internal energy is given by 


S zz {k,oj) 


1 V r —e iult 
2it 

X (&qWV-fcW b g'(°)V+fc(0)) • (5.7) 


2. Magnetic moment and magnetic susceptibility 
The magnetic moment per site is given by 


E=(H)=J2 ~& nk + £ mf, (5.13) 

k 


where the infinite on-site interaction does not contribute 
because its expectation value is zero if the hard-core con¬ 
straint is fulfilled. The specific heat at constant volume 
is obtained by taking the temperature derivative of the 
internal energy, 


C = 


dE 

dT' 


(5.14) 


m = (S z ) 



We note that £mf depends on temperature and has to 
(5.8) p> e taken into account for calculating the specific heat. 
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FIG. 8. (Color online) Contour plots of the spectral function A(k, u>) of the hard-core bosons at k y = 0 for temperatures of 1 K 
and 4 K in a magnetic field B = 9 T corresponding to an energy gap A = 0.19J. The white dashed line is the bare excitation 
energy given by Eq. (2.111. 


fc=0,B = 9T k = Q, B = 9 T 




FIG. 9. (Color online) Spectral function A(k, uo) at k = 0 and at k = Q (at the minimum of the dispersion) in a magnetic field 
B = 9 T corresponding to an energy gap A sa 0.19J. 



FIG. 10. (Color online) Numerical results for the magnetic 
moment m at different magnetic fields between 9 T and 


11.5 T. 


Fig. 13 In Fig. 14 we again compare our numeri¬ 


cal results from the self-consistent ladder approximation 
with the low-temperature Hartree-Fock approximation 
without hard-core interaction and with a simple spin 
mean-field theory described in Appendix C. The mag¬ 
netic contribution to the specific heat of CS 2 CUCI 4 has 
been measured experimentallyP^h^The more recent data 
published in Ref. |T 6 ] differ slightly from Ref. [15] for 
B = 11.5 T and are in better agreement with the ex¬ 
pected size of the energy gap at that field strength. 
Therefore, in Fig. [15] we compare our results with the 
experimental data from Ref. he where we find that 
our theory captures the experimentally observed behav¬ 
ior both qualitatively and quantitatively. At low tem¬ 
peratures, the slope in the logarithmic plot of CT versus 
1/T in Fig. 15 is given by —A, which follows directly 
from Eq. (5.15 I. 


The asymptotic low-temperature behavior of the specific 
heat is given by 

C oc TT e -A/r . (5.15) 

The numerical results for the specific heat at different 
magnetic fields above the saturation field are shown in 


VI. SUMMARY AND CONCLUSIONS 

We have mapped the spin-1/2 Heisenberg model de¬ 
scribing CS 2 CUCI 4 to a model of hard-core bosons where 
the hard-core constraint has been taken into account by 
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T( K) 


FIG. 11. (Color online) Numerical results for the magnetic 
susceptibility \ a f different magnetic fields between 9 T and 
11.5 T. The thin solid lines are low-temperature results from 
the Hartree-Fock approximation without hard-core interac¬ 
tion, which allows us to get results in the low-temperature 
regime where the self-consistent ladder approximation cannot 
be used due to the limited frequency resolution. 


B =9 T 



T( K) 


FIG. 12. (Color online) Comparison of the results for the 
magnetic susceptibility x from the self-consistent ladder ap¬ 
proximation, Hartree-Fock approximation without hard-core 
interaction, and spin mean-field theory for a magnetic field 
B = 9T corresponding to an energy gap A « 0.19J. 


an infinite on-site repulsion. Since we have only con¬ 
sidered magnetic fields B > B c (along the a axis per¬ 
pendicular to the lattice plane), we had to deal with 
gapped hard-core bosons. Due to the energy gap, the 
hard-core interaction can be taken into account using 
the self-consistent ladder approximation 32 and the re¬ 
maining interactions can be treated within the self- 
consistent Hartree-Fock approximation. Before applying 
this method to CS2C11CI4, we have investigated for the ex¬ 
actly solvable one-dimensional XY model how the ladder 
approximation breaks down in the vicinity of the critical 
field B c , finding that the ladder approximation for finite 
energy gaps A works well both at low and high temper¬ 
atures and the deviations, maximal at T ss A, decrease 
with rising energy gap A. We have then calculated the 
spectral function of the hard-core bosons for CS2C11CI4 
from which we have obtained the magnetic susceptibil¬ 



T(K) 


FIG. 13. (Color online) Numerical results for the specific heat 
C at different magnetic fields between 9 T and 11.5 T. The 
thin solid lines are low-temperature results from the Hartree- 
Fock approximation without hard-core interaction, which al¬ 
lows us to get results in the low-temperature regime where 
the self-consistent ladder approximation cannot be used due 
to the limited frequency resolution. 


B = 9 T 



T( K) 

FIG. 14. (Color online) Comparison of the results for the 
specific heat C from the self-consistent ladder approxima¬ 
tion, Hartree-Fock approximation without hard-core interac¬ 
tion, and spin mean-field theory for a magnetic field B — 9 T 
corresponding to an energy gap A ~ 0.19J. 


ity and the specific heat. The calculated specific heat is 
in good agreement with the available experimental data. 
We conclude that the self-consistent ladder approxima¬ 
tion in combination with a self-consistent Hartree-Fock 
decoupling of the non-hard-core interactions gives an ac¬ 
curate description of the physical properties of gapped 
hard-core bosons in two dimensions at finite tempera¬ 
tures. An extension to three dimensions is straightfor¬ 
ward and would only increase the numerical effort due 
to an increasing number of lattice sites. Our meth¬ 
ods can also be directly applied to the material class 
Cs2Cu(Cl4_ a; Br a ,), where chlorine is partially substituted 
by bromine which changes the strength of the exchange 
couplings and the ratio j'/jEHSD \\Tffle in our work 
we started from a spin-1/2 model which we mapped to 
hard-core bosons, our theoretical approach is applicable 
whenever the elementary excitations can be described by 
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1/T (K -1 ) 

FIG. 15. (Color online) Comparison of our numerical re¬ 
sults for the specific heat (solid lines) with experimental data 
from Ref. [16] (symbols) for different magnetic fields between 
9 T and 11.5 T. The dashed lines are low-temperature 
results from the Hartree-Fock approximation without hard¬ 
core interaction, which allows us to get results in the low- 
temperature regime where the self-consistent ladder approx¬ 
imation cannot be used due to the limited frequency resolu¬ 
tion. 


gapped hard-core bosons; some examples are discussed in 
Ref. [32]. In the case of S > 1/2, a mapping to hard-core 
bosons is not known, but mapping the spins to canonical 
bosons is possible by using, for example, the Holstein- 
Primakoff t r ans format ion where the constraint on the 
boson occupation number per site, hi < 2 S, is difficult to 
take into account analytically. Therefore, spin-wave the¬ 
ories based on such a mapping to canonical bosons are 
only valid for low boson densities at low temperatures, 
where the constraint is not important and the usual ex¬ 
pansion in hi/S is justified. 
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APPENDIX A: NUMERICAL DETAILS 

In this appendix, we give more details on the numerical 
solution of the self-consistency equations for the spectral 
function A(k,ui). To find a self-consistent solution for 
A{k 1 uj ), we have to start from an initial spectral func¬ 
tion A init (k,u;). If we would just have the standard non¬ 
self-consistent ladder approximation, we would replace 
A{k , uj) by the non-interacting spectral function Ao(k , w) 
and then directly calculate the spectral function A(fc,w). 
It is therefore sensible to use the non-interacting spectral 
function as the initial spectral function, which is given 

by 

A 0 (k, w) = —ImG'o(fc, w + i0 + ) = S(uj - £*,). (A.l) 

7T 


We note here that at T = 0 the spectral function is not 
affected by interactions, A(fc,u;) = Ao(fc,w), because the 
ground state is the vacuum without any bosons due to the 
energy gap. For the numerical calculation, we replace the 
delta function by a box function of finite width rj (e.g., 
= 0.1 J) centered at w = £*,. This is fine as long as 
£*. > 0 which is the case for magnetic fields B > B c . For 
S,k < 0, we have to take the sign of the spectral function 
into account 

sgn (A(fc, oj)) = sgnw. (A.2) 


Therefore, a positive delta peak is not permitted for neg¬ 
ative frequencies and the non-interacting spectral func¬ 
tion cannot be used for values of k with < 0. In our 
calculations, we instead place a step function at a small 
positive frequency when < 0. This allows us to find a 
self-consistent solution even for B < B c . 

Having chosen an initial spectral function, the next 
step is to calculate p(p, w) via Eq. ( 3.36[ ), which is a 
multi-dimensional convolution that can be calculated 
with the fast Fourier transform method, e.g., using the 
FFTW library.^ Then /(p, ui) can be obtained from 
Eq. (3.401, where the principal value integral can also 
be evaluated as a convo lution.^ Next, the calcula tion of 
p^(k,u)) via Eq. (3.48cI and ReY> R (k,u)) via Eq. (3.48a) 
also involves convolutions. While the values of k are 
naturally discretized for a finite lattice, as discussed in 
Appendix B, the real frequencies to have to be artificially 
discretized, leading to a limited frequency resolution, and 
a frequency cutoff has to be introduced (e.g., |w| < 20 J). 
When using the fast Fourier transform method to eval¬ 
uate the convolutions, this treats the functions as pe¬ 
riodic both in momentum and frequency space, leading 
to a wrap-around effect in the frequency dependence of 
the calculated functions. This wrap-around error can be 
dealt with by setting the spectral function A(k,uj) to 
zero for frequencies larger than a certain cutoff (e.g., for 
|cc| > 10 J). In our calculations, we typically used lattice 
sizes up to 4096 sites and up to 131072 frequency points. 

To achieve convergence, a simple mixing update pro¬ 
cedure has to be used, where the updated spectral func¬ 
tion and Hartree-Fock parameters are set to be a mix¬ 
ture of the previous iteration and the new values from 
the self-consistency equations. In our case, a mixing of 
50% worked well. We note that in the case without the 
self-consistent Hartree-Fock decoupling (e.g., for an XY 
model), mixing is not necessary to achieve convergence. 
The converged numerical result should (up to a small 
numerical error) fulfill the sum rule 32 


duj A(k, ui) = 1 — 2 n. 


(A.3) 


APPENDIX B: BRILLOUIN ZONE 
DISCRETIZATION 

The use of fast Fourier transform methods is based on 
the periodicity of the transformed functions. Therefore, 
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the Brillouin zone should not be arbitrarily discretized 
because that would in most cases destroy the periodicity. 
Still, there is an infinite number of possible parameteri- 
zations of the Brillouin zone. In our work, we have used 
two parameterizations which we will present here. The 
first parameterization starts from the lattice basis 


b c 

a 1 = bx, a 2 = --x+ -y, (B.l) 

with the corresponding reciprocal basis 


where the periodic boundary conditions dictate that 


k[ = 

N{’ 

l[ £ {0,.. 


-1}, 

(B.9a) 

k' 2 = 

1'2 

N 2 ’ 

1' 2 £ {0,.. 


-!}■ 

(B.9b) 

The total number 

of lattice sites is N = 

= N[N' 2 . 

To obtain 

a uniform mesh,®l we 

have to choose 





K = n[. 



(B.10) 


APPENDIX C: SPIN MEAN-FIELD 
APPROXIMATION 


bi 


2lt „ 27T „ 47T „ 

-x-\ -y, b 2 = —y. 

be c 


(B.2) 


The lattice momentum vectors can then be expanded in 
terms of the reciprocal basis, 


We expect that at high temperatures T> J, the spins 
decouple and it is sufficient to describe the spin-spin in¬ 
teractions on a mean-field level where the effects of the 
interactions are approximated by an effective magnetic 
field. To derive this mean-field description, we start from 
the Hamiltonian (2.1), 


k = kibi + k 2 b 2l (B.3) 

where the periodic boundary conditions dictate that 

= Zi e {0,...,^! - 1}, (B.4a) 

k2 = ^~, h £ {0, ...,N 2 - 1} • (B.4b) 

iv 2 

The total number of lattice sites is N = Ni N 2 . To obtain 
a uniform mesh,® we have to choose 

N 2 = 27Vi. (B.5) 

The second (primed) parameterization starts from the 
lattice basis 


n - ! £ [JijSi • .v, + d,j ■ (Si x Sj)] - h^s*. 

ij i 

(C.l) 

First, we note that only the ^-component of the expec¬ 
tation values of the spin operators does not vanish, 

(Si) =mz , m= (S?). (C.2) 

Expanding up to linear order in fluctuations from this 
expectation value, we find 

2 

■Utt-NJo^-hesYtS?, (C.3) 

i 

where the effective magnetic field is given by 

h e s = h — J 0 m, (C.4) 

with 


°-i = 2 X ~2 y ' a2 = 2 X+ 2 Vl 


with the corresponding reciprocal basis 


(B.6) 


J 0 = 2J + 4J'. 


(C.5) 


The magnetic moment m is obtained by solving the self- 
consistency equation 


,, 27t „ 27r . 27 t „ 27t „ 

b l = —x^— y, b 2 =—x+—y. (B.7) 

be be 


The lattice momentum vectors can again be expanded in 
terms of the reciprocal basis, 


k = k ;t>; + k' 2 b' 2 , 


(B.8) 


m = - tanh ( '-h eS ) , 


(C.6) 


and the energy in this mean-field approximation is a sim¬ 
ple function of magnetic field and magnetic moment, 


E = N ( - J 0 m 2 — mh ) . 


(C.7) 
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